package org.hawksee.core.tools;

public class Calculator
{
  public static void invertGauss(double[][] matrix)
  {
    final int n = matrix.length;
    final double[][] inverse = new double[n][n];
    for (int i = 0; (i < n); i++) {
      double max = matrix[i][0];
      double absMax = Math.abs(max);
      for (int j = 0; (j < n); j++) {
        inverse[i][j] = 0.0;
        if (absMax < Math.abs(matrix[i][j])) {
          max = matrix[i][j];
          absMax = Math.abs(max);
        }
      }
      inverse[i][i] = 1.0 / max;
      for (int j = 0; (j < n); j++) {
        matrix[i][j] /= max;
      }
    }
    for (int j = 0; (j < n); j++) {
      double max = matrix[j][j];
      double absMax = Math.abs(max);
      int k = j;
      for (int i = j + 1; (i < n); i++) {
        if (absMax < Math.abs(matrix[i][j])) {
          max = matrix[i][j];
          absMax = Math.abs(max);
          k = i;
        }
      }
      if (k != j) {
        final double[] partialLine = new double[n - j];
        final double[] fullLine = new double[n];
        System.arraycopy(matrix[j], j, partialLine, 0, n - j);
        System.arraycopy(matrix[k], j, matrix[j], j, n - j);
        System.arraycopy(partialLine, 0, matrix[k], j, n - j);
        System.arraycopy(inverse[j], 0, fullLine, 0, n);
        System.arraycopy(inverse[k], 0, inverse[j], 0, n);
        System.arraycopy(fullLine, 0, inverse[k], 0, n);
      }
      for (k = 0; (k <= j); k++) {
        inverse[j][k] /= max;
      }
      for (k = j + 1; (k < n); k++) {
        matrix[j][k] /= max;
        inverse[j][k] /= max;
      }
      for (int i = j + 1; (i < n); i++) {
        for (k = 0; (k <= j); k++) {
          inverse[i][k] -= matrix[i][j] * inverse[j][k];
        }
        for (k = j + 1; (k < n); k++) {
          matrix[i][k] -= matrix[i][j] * matrix[j][k];
          inverse[i][k] -= matrix[i][j] * inverse[j][k];
        }
      }
    }
    for (int j = n - 1; (1 <= j); j--) {
      for (int i = j - 1; (0 <= i); i--) {
        for (int k = 0; (k <= j); k++) {
          inverse[i][k] -= matrix[i][j] * inverse[j][k];
        }
        for (int k = j + 1; (k < n); k++) {
          matrix[i][k] -= matrix[i][j] * matrix[j][k];
          inverse[i][k] -= matrix[i][j] * inverse[j][k];
        }
      }
    }
    for (int i = 0; (i < n); i++) {
      System.arraycopy(inverse[i], 0, matrix[i], 0, n);
    }
  }
}
